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We study the transition time distribution for a particle moving between two wells of a multidi- 
mensional potential in the low-noise limit of overdamped Langevin dynamics. Possible transition 
paths are restricted to a thin tube surrounding the most probable trajectory. We demonstrate that 
finding the transition time distribution reduces to a one-dimensional problem. The resulting tran- 
sition time distribution has a universal and compact form. We suggest that transition barriers can 
be estimated from a single-temperature experiment if both the life times and the transition times 
are measured. 



I. INTRODUCTION 

Dynamical processes in condensed phase often involve 
transitions between metastable states. Although their 
timescales in different systems may differ by several or- 
ders of magnitude, such processes can often be described 
within a common framework. 

A general approach to molecular systems in condensed 
phase is to describe it in terms of an effective potential en- 
ergy represented by a function of the system coordinates. 
The environment, represented by the degrees of freedom 
we are not interested in and cannot control, provides vis- 
cosity, as well as random forces (noise) that make the 
observable dynamics stochastic [lj. If the case of weak 
noise the system spends most of time near the effective 
potential minima and rarely makes transitions between 
potential wells. The above scenario forms a foundation 
of the transition state theory for chemical reaction rates 
dominated by thermoactivated processes. Bulk experi- 
ments, where transition rates are easily observable, mo- 
tivated the development of the transition (reaction) rate 
theory for a variety of situations 0] • The transition rates 
can be obtained directly by observing a chemical reac- 
tion; their ratios can be retrieved from the equilibrium 
concentrations. The probability distribution functions 
(p.d.f.) of these life (dwell) times, i.e., the times the sys- 
tem spends in potential wells between the transitions, are 
represented by the exponential functions, and, therefore, 
are fully determined by the transition rates. 

The bulk experiments, interpreted in terms of the tran- 
sition rate theory, provide some information on the sys- 
tem and the environment. However, the vast majority of 
their properties remain unexplored due to the ensemble 
averaging. The latter, being automatically performed in 
the bulk experiments, hides a lot of dynamical informa- 
tion, e.g., the details of the transition paths 0,0, For 
example, within the low-temperature transition rate the- 
ory, the life times are exponentially long, whereas the 
transitions themselves are treated as almost instanta- 
neous, with their durations being irrelevant. 
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Development of experimental techniques with im- 
proved spatial and temporal resolution (in particular, 
single-molecule experiments) allowed the stochastic phe- 
nomena to be studied with the ensemble averaging 
avoided, thus creating to a new area in statistical physics. 
Theoretical effort has been mostly focused on the fluctu- 
ation theorems and related problems (see, for example, 

Refs. i, 0, S H EE EL Ei El 0, 0, EE Eil)- 

In experiments with single molecules and well- 
controlled systems from nano- to macro-scales [H, 0, 
0, SH 0, SI 0, HE Ull one can observe stochas- 
tic trajectories that drive the system far from equilib- 
rium, even if on longer time or length scales the sys- 
tems exhibit averaged time-independent (stationary) be- 
havior. Detailed information on such processes, e.g., the 
hydrogen-bond rearrangement, can be also obtained via 
the bulk ensemble-averaged measurements, if the proper 
advanced spectroscopic tools, such as 2D infrared spec- 
troscopy [23, SI SI HE SH Hi, SI, are applied. 

One can study different statistical characteristics of 
transition trajectories [H, Si, SI, SI, S3, SI, SI, S3] • In 
this paper we focus on the distribution of the transition 
durations, which is a first passage time (FPT) problem. 
For a given transition trajectory (also referred to as a 
switching path) the transition time is defined is the in- 
terval between the last moment the trajectory leaves a 
neighborhood of one metastable point and the first in- 
stance of time when the trajectory enters a neighbor- 
hood of the other metastable point. A review of work 
on transition-event durations as well as similar problems 
can be found in Ref. [40|- 

Most existing theoretical models that describe multi- 
dimensional transition processes (e.g., in Refs. SE SI, 

0,0, a si) are one-dimensional. Naturally, a realis- 
tic model of a complex system should include more than 
one degree of freedom, and a theoretical background to 
justify such a reduction (identification of a ID reaction 
coordinate) is still being developed. In this manuscript 
we present a systematic study of the transition time 
probability distribution function (p.d.f.) for a multidi- 
mensional potential. In our model the driving force in 
conservative, i.e., represented by the gradient of a po- 
tential function, with the corresponding potential being 
smooth and not having many energy and length scales. 
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Since at low temperatures (or, equivalently, in the weak- 
noise limit) the transitions are rare, one can focus on 
a single transition between two potential wells through 
a saddle point (which may be viewed as a transition 
state). Once the single transition is analyzed, the the- 
ory can be extended to complex realistic systems with 
many metastable states (e.g., hydrogen-bond networks). 

In this manuscript we show that in the weak-noise limit 
the transition time p.d.f. has an almost universal ID 
form. In the above limit the transition time distribution 
is determined by a small number of parameters, related 
to the potential and diffusion tensor. This is due to the 
slowing down of dominating transition trajectories when 
they pass the saddle point. At this point it would be 
worth noting that the ID nature of the transitions be- 
tween the metastable states originates from the effect 
of the transition paths being statistically restricted to 
the narrow tubes that surround the most probable path 
[45[ . It is a general phenomenon that is characteristic 
of low-temperature non-equilibrium and relaxation phe- 
nomena. We have recently demonstrated [46j that in the 
weak-noise limit the statistics of stationary topological 
currents in the infinite-dimensional field theory can be 
found from an effective quasi-lD Markov chain model. 

The details of the transition processes can be especially 
important if their durations, although being shorter than 
the life times (so that the metastable states can be distin- 
guished), are not negligibly short. This is the case of the 
activation barriers being higher than fcgT, yet not ex- 
ceeding several ksT, with ks being the Boltzmann con- 
stant. Such a situation is typical for room temperatures 
and hydrogen bonds, which can be found in various sys- 
tems from biomolecules to bulk and surface water. 

In many cases, e.g. for biomolecules, the properties of 
the thermodynamic (effective) potential and the environ- 
ment are strongly temperature-dependent due to vicinity 
of the melting phase transitions. This substantially com- 
plicates the extraction of the relevant parameters, based 
on measuring the lifetimes at different temperatures. As 
discussed in Sec. [Vj the simple universal form of the tran- 
sition time p.d.f. allows these properties to be retrieved 
from single-temperature experiments, provided both the 
lifetimes and the transition times are measured. 

The manuscript is organized as follows. In Sec. |TT] we 
introduce our basic model of overdamped stochastic dy- 
namics with white Gaussian noise. In Sec. IIIII we define 
the transition duration and calculate the transition time 
p.d.f. Sec. IIVI presents a comparison of the analytical 
and numerical results for moderately weak noise. In Ap- 
pendix we derive an intuitive low-noise relation between 
the transition time p.d.f. and the conditional probabil- 
ity density obtained from the unrestricted Fokker-Planck 
equation. 



II. TRANSITION PATHS 

We consider a system with an m-dimensional config- 
uration space, described by a set r] = (rj 1 , . . . , rf L ) of 
coordinates. The Langevin equation 

il i {T)=F*(rj)+e{r,,T) (1) 

determines stochastic dynamics of the system in the over- 
damped regime. Here F % denotes the deterministic (ad- 
vection) component of the velocity, linearly related to the 
driving force Fj, 

/• ' - .'/"/ •, , (2) 

via the mobility tensor g tJ {r]). The deterministic dy- 
namics is assumed to be conservative, i.e., Fj — —djV, 
with V(r]) being the corresponding potential. In this 
manuscript we focus on a situation when the potential 
V{rj) is bistable, i.e., it has two local minima (basins) 
V{Vi) = ^i an d y{V2) = V%. The ways of extension of 
our approach to a more general situation with multiple 
basins is briefly outlined in Sec. [V] 

The mobility tensor can be viewed as a Riemann metric 
in the configuration space. Due to the Einstein relation 
(fluctuation-dissipation theorem), which reflects the fact 
that the bath is at equilibrium with the temperature T, 
the same tensor g u weighted with a factor k = ksT 
characterizes the correlations for white Gaussian noise: 

(e(V,r 2 )e(v,rx)) = KgV( V )6(T 2 - n ). (3) 

The dynamics represented by Eqs. (JTJ) and allows for 
a path integral description with the Onsager-Machlup 
[47| action (hereafter we imply summation over the re- 
peated indices) 

S( V ) = (1/2) [ drg ik {ff - F*(r,)) (r, k - F k ( V )) , (4) 
Jo 

where gij is defined by gijg^ k = S k . The transition prob- 
ability K(x" , x'; t), i.e. the probability for the system to 
move from point x' in the configuration space to point 
x" over time t is given by a path integral 

K(x", x'; t) = J Vr]e- S ^I K , (5) 

where integration goes over all trajectories with 77(0) = 
x' and rj(f) = x" . The action in the path-integral repre- 
sentation [Eq. |5])] is understood as a proper discrete-time 
action, whose continuous limit is given by Eq. Q [48]. 

In the weak-noise k — > limit the expression for 
K(x" ,x';t) [Eq. (J3J)] suggests that a transition is domi- 
nated by the trajectory that minimizes the action S(rj). 
It is well-known that integration in Eq. goes over 
continuous, rather than smooth trajectories [49(, which 
can be actually considered as a reason why the action in 
Eq. defined for smooth trajectories, requires regu- 
larization (i.e., a choice of a proper discrete-time form). 
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However, in the limit of a fine discretization (small time 
step) the dominant trajectory becomes smooth, and, 
therefore, can be found from the Onsager-Machlup action 
@ by solving the corresponding Euler-Lagrange equa- 
tion. 

There are two special types of trajectories that provide 
local minima of the Onsager-Machlup action. Its abso- 
lute minima (S = 0) correspond to the trajectories that 
satisfy if = F l (r](t)). However, for such downhill trajec- 
tories the value of the potential decreases monotonically 
with time, whereas a transition between two basins in- 
volves the potential increase at the first stage of the pro- 
cess. In the conservative Fj — —djV case under consid- 
eration an uphill trajectory that satisfies rf = —F' l (r](t)) 
also satisfies the Euler-Lagrange equation. Note that 
an uphill trajectory is nothing else that a time-reversed 
counterpart of a downhill trajectory. 

The scenario of a transition between two basins has 
been thoroughly studied in the context of the reaction 
rate theory [lj. The system follows the uphill trajectory 
to reach the transition state (a critical point of the po- 
tential with a single unstable mode) and then reaches 
the other basin via the downhill trajectory. The uphill 
and downhill parts coming to rest at the saddle point 
constitute a special trajectory (hereafter referred to as 
the reference trajectory) that minimizes the Onsager- 
Machlup action. The downhill part does not need any 
noise, whereas for the uphill part certain concerted noise 
is necessary to overcome the deterministic force. The two 
parts of the reference trajectory satisfy the equations 

f,*=±F*( V (t)). (6) 

The reference trajectory is a special case of the most 
probable escape path (MPEP) or the most probable 
switching path (MPSP) that corresponds to infinite tran- 
sition time. We assume that the reference trajectory 
connecting two potential minima is unique. When the 
transition time is finite, the most probable trajectory de- 
viates from its reference counterpart. In particular, in 
the limit of very short times it becomes a geodesic line, 
and is independent of the potential. 

Within the path integral approach we need to consider 
all possible transition trajectories. In the weak-noise case 
the situation is much simpler, since the deviations from 
the reference trajectory are suppressed. This can be 
used, for example, for generalizing ID results [5(J HH 
to the multidimensional case, by linearizing the theory 
around the reference trajectory. Thus, we have a lin- 
ear but still multidimensional time-dependent problem. 
Therefore, the deviations of the MPSP from the reference 
trajectory due to finite transition times, as well as Gaus- 
sian fluctuations around the MPSP may not be treated 
explicitly within the analytical framework. 

To avoid the aforementioned difficulties we invoke 
a "Hamiltonian" , rather than "Lagrangian" approach: 
Instead of following individual stochastic trajectories, 
which is the natural basis for the path-integral approach, 
we reformulate the problem in terms of the probability 



density evolution. The probability dP to find a system 
in a certain volume of the configuration space is obvi- 
ously an integration measure, and does not depend on 
the particular choice of coordinates we use to parame- 
terize the configuration space. Therefore, the probability 
density is a scalar if defined as p = dP/ dp with respect to 
the invariant volume element dp — ^fgdx 1 . . . dx m , where 
g = det gik. With this invariant definition of the proba- 
bility density, the Fokker-Planck operator, which governs 
its evolution, is naturally invariant with respect to coor- 
dinate transformations. 

The Langevin equation ([T]) corresponds to the follow- 
ing Fokker-Planck (FP) equation: 

8tP = £p, (?) 
where the FP operator has a covariant form [53| 
Cp = {n/2)V 2 p - div (Fp) = 

(«/2)(i/V5)ft (Vgg* k dkp) + (V^V) ■ (8) 

III. TRANSITION TIME DISTRIBUTION IN 
THE WEAK-NOISE LIMIT 

A problem of the transition duration and its distri- 
bution can be formulated for any noise strength. For 
strong noise the framework of overdamped dynamics with 
white Gaussian noise, introduced above, may be inade- 
quate, and a model with more parameters might be nec- 
essary. Therefore, we restrict our treatment to the low- 
noise limit. 

In this limit a conservative system has the following 
properties [3, H3, [5f| • First, the quasiequilibrium den- 
sity distributions are strongly localized in the potential 
wells, and the inter-well relaxation requires exponentially 
long times. The shapes of these distributions are deter- 
mined by a small number of parameters that characterize 
the potential in the wells. The relaxation rates, in ad- 
dition, depend on the mobility and the potential at the 
transition states represented by the saddle points. Sec- 
ond, the trajectories that lead from one potential well to 
another are also determined by a few parameters of the 
deterministic drift, as we will see below. This universality 
that takes place in the low-noise limit makes the model 
especially simple and attractive, even though it can be 
imprecise when applied to stronger noise. 

The analysis of the Onsager-Machlup action or the FP 
operator shows that a characteristic transverse deviation 
of the trajectory from the most probable trajectory is oc 
y/n. For weak noise, most observed transition trajectories 
are close to the reference trajectory, since they typically 
correspond to long transition times, as follows from the 
results of Subsection IIII Dl Thus, all probable transition 
trajectories can be considered to be enclosed in a thin 
tube of the width oc \JH. The tube extends from one 
potential well to the other and is represented by a tubular 
neighborhood of the reference trajectory. The concept of 
a tube, as a region outside of the potential minima where 
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the transition density is localized, appeared, for example, 

in Refs. HE H m EE M, HE m 



A. Coordinates in the transition tube 

Once the reference trajectory [the MPSP correspond- 
ing to the infinite transition time defined in Eqs. ([6])] 
in the m-dimensional configuration space is known, we 
can introduce a coordinate system in the tube that al- 
lows all possible transition trajectories to be described 
in a convenient and natural way. The coordinates are 
represented by the distance x along the reference trajec- 
tory and a set of (m— 1) mutually orthogonal transverse 
coordinates We set x — at the saddle point and 
denote the longitudinal coordinates of the first and sec- 
ond minima by y\ < and j/2 > 0, respectively. The 
transverse coordinates can be selected in the following 
way. For a point on the reference trajectory located at 
x and a tangent vector £ (i.e., an infinitesimal shift of 
a configuration) at this point, which is normal to the 
reference trajectory, we consider a length |£| segment of 
a geodesic line (recall that the mobility tensor defines a 
Riemann metric in the configuration space) that starts at 
the point x on the reference trajectory in the direction 
of £. The end point of the segment is denoted by (x,C). 
For small enough |£| (actually small compared to the tra- 
jectory curvature radius) the sets (x, £) unambiguously 
identify the configurations in the tubular neighborhood 
of the reference trajectory. At this point a global coor- 
dinate system in the tube can be built by choosing the 
basis sets e\{x), . . . , e m ^i(x) in the orthogonal spaces to 
the reference trajectory with a smooth dependence on x. 
Naturally a coordinate set (x, C 1 , ■ ■ ■ , C™ 1 ) represents 
the point (x, C) with C = Y,a=i C a e a (x). 

Since we are dealing with a weak-noise (low- 
temperature) and long-time situation, so that the devia- 
tions from the reference trajectory are small, we can ex- 
pand all relevant quantities in the deviations £ from the 
reference trajectory. Thus, the potential in the vicinity 
of the tube is given by 

V( V ) = Vo(x) + (\/2)W aP {x)CC P • (9) 

The terms linear in £ do not appear in the expansion, 
since the potential gradient is directed along the reference 
trajectory. In the low-noise limit the expansion to second 
order in £ is sufficient, provided W is non-degenerate, i.e., 
W is positively defined along the transition trajectory. 
Therefore, we will need the potential only the form given 
by Eq. ([9|), and in particular use the longitudinal force 
on the reference trajectory 

F Q (x) = -d x V (x) , (10) 

which has zeros in the stationary points yi, 0, and yi. 




FIG. 1: Most transition trajectories between neighborhoods 
of the stable points are within a transition tube. 



B. Definition of transitions 

We are focused on a situation when the particle makes 
a transition from one potential well to another via a single 
"channel" that passes through the potential saddle point. 
We assume the temperature (noise) to be small compared 
to the height of the potential hill. Thus, a particle located 
in a potential well has enough time to reach an almost 
equilibrium distribution before a rare strong-noise fluctu- 
ation, represented by a series of concerted random force 
kicks, pushes it over the potential saddle point to the 
other potential well. A transition is naturally defined as 
an event when a particle leaves a neighborhood U\ of the 
initial basin, located at x = y\, and enters a neighbor- 
hood U2 of the other minimum located at x = y% . Thus, 
the transition time is the first passage time [52| from the 
initial boundary dU\ to the final boundary dXJi without 
revisiting the initial boundary. It is convenient and nat- 
ural to choose the regions U\ and U2 in such a way that 
the probability to find the particle outside of them is 
small. On the other hand, since we study the transitions 
between metastable states, the regions Ui and U2 should 
be small compared to the typical length scale of the prob- 
lem. These two conditions are compatible, provided the 
noise strength n is small, compared to the barrier height. 
The subsequent treatment of the problem can be sim- 
plified without losing the essential properties of its solu- 
tion, if we assume the faces of the tube at x = Xi and 
x — X2 to be parts of the boundaries dU\ and 0U2, re- 
spectively (see Fig. [IJ . Therefore, the explicit conditions 
for x\ and X2 read: y 1 —k/ d x F^(yi) <§; (x\ —yx) <C 1 and 
\/—n/d x F (y2) <C (2/2 — £2) *C 1, where the coordinate is 
rescaled, so that the typical length scale of the potential 
can be denoted by 1. 

The transition time p.d.f. can be represented as a path- 
integral form, 

v(t)edu 2 

F[t) = J Pr/e~ s ^ , (11) 

with the integral measure including only those trajec- 
tories r}(r) that satisfy the initial and final conditions, 
specified in the integral and do not visit U\ or U2 at 
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r £ (0, i). As already mentioned in Sec. [II] the path in- 
tegral can be defined once a proper time discretization is 
introduced. If time discretization involves n intermediate 
time points between and t, the regularized integration 
is represented by n integrals over rj £ M\(UiUU2), an in- 
tegral over dU% with a proper distribution of the starting 
points, and an integral over dXJi. The presented above 
discretized form of the path integral has a direct interpre- 
tation in terms of the experimental observables. Namely, 
in an experimental set-up one can monitor the particle's 
position at certain discrete times and register the cross- 
ing of the boundary if the particle appears on its other 
side after a time step. 

We focus on the transition from one basin to the other; 
in the low-noise case they occur by passing the relevant 
saddle point. Since the noise is weak compared to the 
deterministic force on the tube faces at x\ and X2, the 
probability of the excursions into U\ or U2 is low for a 
trajectory that starts at dU\ and ends at dU^ after a 
relatively long time. Therefore, we can approximate the 
expression (fTT|) . which is exact but difficult to deal with, 
by the matrix element 

F(t) ex {p'\e ct \ e ) (12) 

with the Fokker- Planck operator C given by Eq. (|SJ . The 
proportionality sign means that the distribution should 
be normalized, since not all trajectories that start at dU\ 
will reach dU% (actually, only an exponentially small frac- 
tion of them will do). However, as explained above, the 
trajectories that do reach dU2 pass through the tube of 
the small width ~ y/H and mainly stay within the inter- 
val (3:1,0:2). In Eq. [)12p \g) represents the distribution 
of the starting points at d\J\ , and (p' | stands for the inte- 
gration over the final surface. The initial distribution has 
an equilibrium form and at x = x\ is contained within 
the face of the tube: 

0(17) = Ze- K ~ lw{v ^®Q 5 ( x _ (13) 

where we made use of the fact that x\ is close to the 
minimum at x — y±. 

In Appendix [X] we derive the approximation repre- 
sented by Eq. (fT2]) in the ID case, which can be gen- 
eralized to a multidimensional situation. 

C. Eigenstates of C 

The matrix element in Eq. (|12|) can be found as a spec- 
tral decomposition in terms of the eigenstates of the FP 
operator, which are the solutions of the eigenvalue prob- 
lem Cp = Xp. The calculation below can be viewed as a 
simplified version of the multidimensional calculation re- 
ported in Ref. 0, adapted to the initial distribution near 
the potential minimum. We skip some finer details (irrel- 
evant in the low-noise limit) of the derivation replacing 
them with physical intuition; these details can be found 
in Refs. QUI 



We seek for the eigenstates in the form 

p{<q) = p (x) exp(- K "V(x)(C ® 0) (14) 

with an x-dependent matrix a. Together with the noise 
strength k, the matrix a determines the density distribu- 
tion in the transverse direction. The resulting eigenvalue 
problem for po(x) is represented by 

( K /2)d 2 xPo - d x ( Po F ) - p Tr(a -W) = X Po . (15) 

There are two kinds of regions along the tube: a do- 
main where the diffusion terms are small compared to 
the deterministic counterparts, hereafter referred to as 
the WKB domain (using the analogy [1^ between the 
Fokker-Planck and Schrodinger operators in the conser- 
vative case under study), and the fluctuation regions lo- 
cated near the potential stationary points. In the weak- 
noise k — ► limit, the size of the fluctuation regions is 
proportional to a/k. Therefore, the WKB domain over- 
laps with the region |a;| <C 1 near the saddle point where 
the force is linear in the deviations (or the potential is 
harmonic, so that the fluctuation regions can be natu- 
rally referred to as the harmonic regions) . In addition to 
Eq. (fT5"|) . the eigenvalue problem Cp = Xp leads to an 
equation for cr(x), which has a simple form in the WKB 
domain: 

2cr 2 — aW — Wa + FqV x <j = . (16) 

Here V x denotes a long (covariant) derivative with re- 
spect to the connection that takes into account the fact 
that the tensor a(x) is defined in different spaces at dif- 
ferent x. Note that the configuration-dependent mobility 
tensor (metric) implicitly enters Eqs. (fTP)) , (fTS")) , and (|T5|) 
through the matrix operations and the longitudinal coor- 
dinate x chosen as a length. In the derivation and results 
below we use the potential in terms of x. 

We will avoid explicit analysis of a outside the WKB 
domain and use its following properties without a rigor- 
ous proof (see, for example, [46j for some details). Close 
to the critical points of the potential we neglect devia- 
tions of a from W at the corresponding fixed points. In 
particular, in the saddle point er(0) = W(0), and on the 
faces of the tube a(x\) w W(yi) and £r(x 2 ) ~ Wfa). 
Thus, the initial distribution (flU]) corresponds to the 
ground state in the transverse direction represented by 
a Gaussian. This excludes the transverse excited states 
from the expansion of e ct g over the eigenstates of C. In 
addition, transverse excited states do not make any con- 
tribution because the corresponding contributions vanish 
upon the integration performed at dU2- 

At this point we would like to emphasize that although, 
as shown above, the higher transverse modes do not con- 
tribute to the transition time p.d.f., the full transition 
probability K(x" , x'; t) involves the complete expansion, 
and the higher transverse modes are not suppressed by 
any small parameter. This is the "Hamiltonian" signa- 
ture of the fact, briefly discussed at the end of Sec. [Ellin 
terms of the "Lagrangian" path-integral language, that 
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although the transition amplitude is dominated by the 
most probable trajectory, it may not be found explicitly 
in analytical terms due to the complexity of the relevant 
multidimensional linear problem. This demonstrates the 
advantages of the "Hamiltonian" approach based on the 
Fokker-Planck equation. This also demonstrates the non- 
trivial nature of the ID-reduction of the transition-time 
distribution problem, as resulting from exact cancela- 
tions related to higher-order transverse modes. 

Next, we notice that the choice of X\ and x 2 described 
in the subsection (|III B[) leads to the same first-order 
equation for a for each longitudinal mode (that corre- 
sponds to the ground-state transverse mode), valid al- 
most everywhere on the interval (xi,X2), as we can see 
from Eq. (Tl6|) . Therefore, <r(x) is the same for all 
eigenfunctions; it is represented by the solution of the 
same equation ([16]) with the same boundary condition 
a{xi) — W(yi). These arguments justify our ansatz rep- 
resented by Eq. (TT4")) . 

The differential operator in Eq. (TTS]) is not self-adjoint 
(Hcrmitian) but can be made self-adjoint by applying the 
transformation 



p (x) = ip (x) exp(-K 1 V Q (x)) . 



(17) 



The transformed eigenfunctions ipo( x ) satisfy the equa- 
tion 

(k/2)8 2 Jj - (/s-72)^ (k8 x F + F 2 + Tr(a - Wj) 
= A^o (18) 

that can be viewed as the Schrodingcr equation with an 
effective potential with deep minima near the fixed points 
(stable, as well as unstable) of V(rj). The eigenfunctions 
responsible for transition dynamics are concentrated at 
the minimum of the effective potential that corresponds 
the saddle point of V(rj). 

To calculate the transition time p.d.f. we expand the 
matrix element (Ti"2"j) of the FP evolution operator in the 
eigenmodes of C as 

(p'\e ct \ ) = Ze K ~ 1( y°^- Vo ^ 



/ (ttk)™" 1 
det¥(^) 



^VoA(xi)VoA(x 2 )e At . (19) 



Note that only the eigenstates localized in the transition 
region provide non-negligible contributions to the expan- 
sion, since only these eigenstates are "produced" by tran- 
sition trajectories that do not visit the neighborhoods of 
the potential minima. 

The next step requires solution of the eigenvalue prob- 
lem, given by Eqs. (fT8|) and (|16| . Since the relevant 
eigenstates are localized at the saddle point, it is suffi- 
cient to consider the harmonic region that contains the 
saddle point x = 0, as well as the two surrounding (and 
overlapping with it) WKB sub-domains. 

In principle, in each WKB sub-domain the eigenfunc- 
tion ^qa contains two "waves" that correspond to two 



linearly independent solutions of the second-order dif- 
ferential equation. One "wave" decays as \x\ increases, 
whereas the other one grows. For the eigenstates un- 
der study (localized near the saddle point) the grow- 
ing "wave" can be neglected immediately This implies 
that functions poA = ipo\ e ~ Vo< " x ^ K change slowly with x. 
Therefore, the shortest way to obtain the solution is to 
neglect the second derivative in Eq. (TT5|) and use Eq. 
(|16p to exclude a(x). The resulting first-order equation 



d x (F Q p 0X /Vdeta) = -Xp ox /Vdeta (20) 
can be easily integrated, which involves one constant: 

dz 



Vdet cr(x) 
Po\(x) = D x — — — — exp -A 
Fa{x) 



Fo(z) 



(21) 



In the harmonic region \x\ <C 1, Eq. (|18p can be rep- 
resented as 

(K/2)d^ - ^k 2 x 2 /{2n) = (A + fc/2)V>o , (22) 

where we approximate a(x) as W(x) = W(0) and denote 

k = d x F (0). (23) 

Since in the WKB domain there is only one "wave" that 
decays as \x\ increases, the eigenfunctions in the har- 
monic region have the form 



ipa 



■)2-2nl 



■kk,(ti — l)! 2 



1/4 



-kx 2 /(2K) 



(24) 



whereas the eigenvalues are given by A = — nk, with in- 
tegers n > 1 and Hermite polynomials H n . The normal- 
ization factor in Eq. (|24[) is determined by the form of 
the solution in the harmonic region where it is localized. 

Now we can match the solutions in the two WKB sub- 
domains with the solution in the harmonic region. We 
will only consider the second WKB sub-domain from x ~ 
^J~k to x — X2; the treatment of the first WKB sub- 
domain is similar. If we choose X2 as the lower limit of 
the integral in Eq. (|2"Tj). the WKB solution with A = — nk 
can be represented as 



Pox{x) - D x k- 1 \/detW(Q)x n ' 1 (K/x 2 ) n , (25) 



K = cxp\ / dx(l/x- k/F (x)) 



(26) 



at x <gC 1. This solution should be matched with the solu- 
tion p x(x), given by Eqs. (fl"7|) and ([24]) . at x » y/n/k in 
the harmonic region. To achieve that, we use the asymp- 
totic form H n (z) = (2z) n of the Hermite polynomial, for 
z > 1. 
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Eventually we obtain the following expression for the 
eigenfunction values at the right boundary: 



1/4 



2fcV"' 1)/2 k{x 2 /KY 
k) ^/(n-l)\F (x 2 ) ■ 

The expressions for po\(x\) are derived in a similar way. 

D. Calculation of the transition time p.d.f. 

We are now in a position to complete the calculation of 
the transition time p.d.f. We first note that its long-time 
asymptotic form can be easily found once we know the 
lowest eigenvalue of — C. The corresponding eigenfunc- 
tion is localized in the vicinity of the saddle point, where 
the potential is harmonic. Therefore, the asymptotic 
form F(t) oc e~ kt is fully determined by the curvature 
k of the potential Vq (the derivative of the longitudinal 
deterministic drift) at the saddle point. However, since 
the expansion over the eigenstates of C is singular in re 
(as shown below), this first term in the asymptotic ex- 
pansion may not be used to determine F(t) in a broader 
region, i.e., around the most probable transition time. 

We complete the derivation by collecting all terms in 
the expansion (flU)) . The time dependence of the transi- 
tion time p.d.f. is given by 



OO 

F(t) ocY,- } (-2fcG/«) n e- fe (" +1 >* 



(28) 



n=0 



where G depends on x\ and x 2 as 



/ X2 \ 

G = fc|xi|x 2 exp J dx(k/\F (x)\ -l/\x\) 



(29) 



We further notice that the series in Eq. (|28|) represents 
a Taylor expansion of an exponential. After normalizing 
the result for a small re, we obtain the transition time 
distribution in a form 



F(t) = (2fcG/re)exp(-fci- (2G/re)e~ fct ) 



(30) 



The form of the obtained expression is the same as the 
one that would be obtained in a ID transition time prob- 
lem. The propagator for the ID FP equation in the 
weak-noise limit has been calculated in Ref. [5l| by us- 
ing a semiclassical approximation for the corresponding 
path integral, yet we are not aware of the simple form 
(|30|) being discussed even in connection with a purely 
one-dimensional transition time problem. 

The result (|30"|) is asymptotically correct for re — > 0. In 
realistic systems, where the details of the transitions can 
be observed, re may not be too small. In the next section 



we present some numerical results for moderately weak 
noise, and demonstrate that our analytical theory still 
provides an adequate picture, at least on the qualitative 
level. 



IV. NUMERICAL RESULTS 

In this section we compare the estimates, provided by 
our approximate theory, with the numerical results for 
the transition time p.d.f. We restrict our comparison 
to ID examples, since the distribution function in the 
weak-noise limit has the same form as in the ID case. If 
our multidimensional potential does not possess any ad- 
ditional length or energy scales, and the potential curva- 
tures in the transverse directions are not atypically small, 
the numerical results presented below are characteristic 
of a general multidimensional case. 

To calculate the transition time p.d.f. numerically we 
use the following scheme [13, H2|- We first solve the 
time-dependent FP equation supplemented by the fol- 
lowing initial and boundary conditions. The initial con- 
dition corresponds to injection of particles close to the 
left boundary at x — X\ + a: 



p{x, 0) = 5{x — x\ — a) 



(31) 



with a small parameter a. The boundary conditions are 
absorbing: 



p(xi,t) = p{x 2 ,t) = 0. 



(32) 



The solution p(x, t) determines the probability flux 
Lo(t) = —(n/2)d x p(x 2 ,t) through the right boundary, 
which is proportional to the transition time p.d.f. Fi- 
nally, one normalizes the flux and applies the a — ► 
limit: 



F(t) = lim 



°fZ°dtu>(t) 



(33) 



Note that the flux through the right boundary vanishes 
when the starting point approaches the left boundary, 
lim a ^ u>(t) = 0. 

The comparison of the transition time p.d.f. obtained 
numerically and from Eq. (|30[) is shown in Fig. [3] for 
an anharmonic potential sketched in Fig. [31 We use the 
corresponding harmonic potential (also shown in Fig. [2]) 
to estimate the source of discrepancies in the transition 
time p.f.s.'s. Since the propagator of the unrestricted FP 
equation (which in ID coincides with the matrix element 
in Eq. [|12[0 is known exactly for the harmonic potential, 
it is easy to check the approximation (|12jl for the transi- 
tion time p.d.f. by the relevant matrix element, as well as 
the subsequent approximation for the the matrix element 
by the spectral expansion (|28|) . Both approximations are 
of the WKB type: they imply weak fluctuations. In fact, 
we observe that the accuracies of the two approximations 
are similar. The final expression (|30|) even turns out to 
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FIG. 2: Two potentials used to illustrate the theory: an 
inverted parabola Vo = x 2 /2 and an anharmonic potential 
Vo = 6.25 sin(i/2.5). Transition boundaries are xi = —5 and 

X2 — 5. 



k, 0.3 
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FIG. 3: Transition time distributions for the anharmonic po- 
tential shown in Fig. [2] from Eq. (|30[) (lines) and numeri- 
cal results (lines with symbols). Distributions are shifted to 
longer times as the noise strength decreases. 



be slightly better than the approximation (TT2")) when it 
comes to such characteristics of the transition time distri- 
bution as the most probable and average transition times. 
A noticeable deviation of the expression ([50)1 from the nu- 
merical results occurs at stronger noise, when the former 
overestimates the probability of both the shortest and the 
longest transition times. At moderately weak noise our 
expression (|30|) correctly approximates the overall shape 
of the distribution function, as well as its parameters. 

Figure [4] shows the first two eigenvalues of the FP op- 
erator £ with the absorbing boundary conditions. The 
lowest eigenvalue of — C determines the asymptotic decay 
rate of the transition time p.d.f. The second eigenvalue, 
in particular, estimates the region where the asymptotic 
regime is valid. The zero-noise limits for the eigenvalues, 
which enter the expansion {25]), are shown in the figure 
with dotted lines. 

The spectrum of the FP operator becomes strongly de- 
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FIG. 4: Two lowest eigenvalues of — C with absorbing bound- 
ary conditions for two potentials shown in Fig. [2] harmonic 
(red squares and circles) and anharmonic (black triangles). 



pendent on the boundaries x\ and xi when they exit the 
WKB sub-domains surrounding the saddle point. This 
is reflected in the overall shape of the transition time 
p.d.f. For instance, if the starting point x\ moves to 
the left beyond the fluctuation region near the poten- 
tial minimum at x = yi, the corresponding "transition 
time" obviously includes exponentially long (if noise is 
weak) residence near the potential minimum; the result- 
ing "transition time" distribution resembles the life time 
distribution with the exponentially small asymptotic de- 
cay rate. 

The shifts of the distributions to longer times, which 
we see in Fig. [3J can be easier observed when analyzing 
the behavior of the most probable and average transition 
times. These characteristics of the p.d.f. are presented 
in Fig. [5] The agreement between the approximation 
([50)) and the numerical results is satisfactory even for 
moderate noise strengths. The quantitative agreement 
naturally improves when the noise becomes weaker com- 
pared to the deterministic force outside the saddle point 
vicinity. 

The result in Eq. (|30|). asymptotically correct in the 
weak-noise limit, is valid for t ^S> k^ 1 . Most importantly, 
the region of validity includes the vicinity of the maxi- 
mum of F(t) at <o = k Ih(2G/k). The average tran- 
sition time is given by (t) = k~ 1 (-fE + ln(2G/re)) with 
-fE ~ 0.577 being the Euler constant. The transition 
times only weakly depend on the potential far from the 
saddle point. As long as k remains small (k <C G), its 
variation results only in the shift of the time that pre- 
serves the shape of the p.d.f. In particular, the variance 
of the transition time (< 2 ) — (t) 2 = n 2 /(6k 2 ) does not 
depend on the barrier height. Our result for the transi- 
tion time p.d.f. is independent of the transition direction, 
which is a manifestation of the time-reversal symmetry, 
which is in place also for the case of more general, not 
necessarily overdamped, stochastic dynamics [39j . 
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FIG. 5: Comparison of characteristic times obtained from 
numerical solution of the FP equation (lines) and from the 
approximate p.d.f. (|30|) (symbols) for two potentials shown in 
Fig. [2] harmonic (dashed line, red squares) and anharmonic 
(solid line, black triangles). Top: the most probable transition 
times; bottom: the average transition times. 



V. DISCUSSION 

We have calculated the distribution of the transition 
time for a particle moving between two stable points 
in a multidimensional potential field, when the noise is 
low compared to the potential hill. We showed that the 
stochastic paths that dominate the transition process are 
restricted to a thin tube. This reduces the calculation of 
the transition time p.d.f. to an almost ID problem. 

Our main result, expressed by in Eq. (I30p . which is 
asymptotically correct in the weak- noise limit, is valid for 
the transition times t ^> fc _1 . We do not consider shorter 
transition times. First of all, they require stronger noise 
intensity, which can invalidate the description of dynam- 
ics by the overdamped Langevin equation with white 
Gaussian noise. Second, even if the white Gaussian noise 
description is valid, the distribution looses the universal 
form at shorter times (which are more difficult to ob- 
serve), which substantially complicates utilization of this 
theory for the interpretation of the experimental data in 
terms of the underlying effective potentials. 

The divergence of the transition times in the limit k — > 
is a result of the competition between the deterministic 



and random forces. In this regard, in the weak-noise limit 
the transition problem we consider is similar to the of 
escape from the saddle point 0, [5{| . The asymptotic 
form of the transition time p.d.f. is also similar to that 
of the escape time distribution, cx e~ kt . In both cases 
its form is determined by the lowest eigenvalue of the FP 
operator corresponding to the eigenstate localized near 
the saddle point. 

The experimental measurements of the transition 
times can be utilized to obtain more information on the 
effective free energy. In the expression for the transi- 
tion rate 7 = C{j3) exp(—(3E((3)) both the prefactor and 
the energy barrier typically depend on the thermal en- 
ergy /3 _1 = ksT. In particular, in the overdamped 
model the prefactor for the transition from the poten- 
tial minimum at y\ through the saddle point at x = 
is given by C = y/kkx det Wtjjij/ det W(0)/(2n), where 
k\ = —d x F (yi), k — d x Fo(0), and W(x) characterize the 
potential V(x, £). The latter is is obtained from the phys- 
ical potential by rescaling the coordinates with the mobil- 
ity tensor; both the physical potential and the mobility 
tensor being temperature-dependent. Therefore, E(f3) 
may not be found by just measuring the rates at differ- 
ent temperatures. In biological molecules, for example, 
the temperatures of interest often lie in a narrow interval 
close to the first-order phase transition, which leads to a 
pronounced temperature dependence of C((3) and E(f3). 
In many cases (such as folding of small RNA fragments) 
a typical lifetime of a metastable state varies from hun- 
dreds milliseconds to several seconds, while the barriers 
are known to be relatively low {e~^ E ~ I0 _1 -r 10 -3 ). 
This suggests that the mobility has atypically low values 
and hence can strongly affect the transition rates. Thus, 
the rate measurements alone cannot distinguish the ef- 
fects of lower mobility in the preexponential factor C ((3) 
from those that originate from higher barriers. In the 
simple overdamped model presented in this paper, if cur- 
vatures of the potential V(f]) in different directions in 
the minima and in the saddle points have similar abso- 
lute values, C can be estimated as equal to this typical 
curvature. The same quantity determines the asymptotic 
decay rate in the transition time p.d.f. Thus, a measure- 
ment of the transition time p.d.f. and the transition rate 
at a single temperature T = l/(fcg/3) allows the barrier 
E((3) to be estimated. Additional information on the po- 
tential landscape and the mob ility can be obtained from 
the intra- well relaxation time [24j . 

In this manuscript we studied in detail the properties 
of transitions between two metastable states. The results 
for the transition time distribution were obtained in the 
limit of vanishing noise. In realistic systems where the 
transitions can be observed, noise cannot be too weak. 
We expect that although under such circumstances our 
results become less precise, they still adequately describe 
the general features of the transition time distribution. 

The developed approach can be extended to a mul- 
tidimensional system with multiple metastable states, 
in the weak noise (low temperature) limit, the result- 
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ing model is reduced to a Markov chain process, de- 
fined on a graph, whose links represent the saddle points 
(transition states), and are described by the transition 
time distributions. An interesting example of such a sys- 
tem would be a hydrogen-bond network [60j whose dy- 
namics can be observed using 2D infrared spectroscopy 
[H IH [H, [H HI M, [H. Although the H-bond rear- 
rangements in water happen on much faster time scales 
compared to the conformational transitions in typical 
biomolecules, the corresponding transition paths can be 
still observed indirectly by studying evolution of the vi- 
brational frequency, available via 2D spectroscopic mea- 
surements H|| . Application of the model described in 
this paper to hydrogen-bond networks will help to con- 
nect the transition mechanisms and transition state prop- 
erties in these complex multidimensional systems to the 
observed spectroscopic signatures. 
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APPENDIX A: REDUCTION OF THE P.D.F. 
PROBLEM TO CALCULATION OF THE 
UNRESTRICTED PROPAGATOR 

In this Appendix we show that if the transition oc- 
curs between two WKB sub-domains, the transition time 
p.d.f. can be approximated by the normalized matrix el- 
ement of the evolution operator e ct , as suggested in Eq. 
(|12p . This property expresses the extremely low prob- 
ability of the deviations from the most probable unre- 
stricted trajectory in the low- noise limit. The deriva- 
tion below is performed for the ID case, since only the 
motion along the reference trajectory is relevant in the 
weak-noise limit for any number of dimensions. The ID 
counterpart for the matrix element, introduced in the 
main text, is the propagator P(x, t) for the unrestricted 
FP equation (i.e., the equation without specific bound- 
ary conditions), whose spectral expansion contains only 
those eigenfunctions that correspond to the transition un- 
der study. 

Consider the solution P(x, t) of the FP equation with 
the absorbing boundary conditions (|32[) and the initial 
condition (|3 1[) . The transition time p.d.f. can be ob- 
tained as the limit of the normalized current through 
the right boundary generated by P(x, t) according to Eq. 



r \t) — — lim — — — r- , 

a "*° L„ I" Vo/A^i + a)d x ifQ^{x 2 ) 



(Al) 



absorbing boundaries, tpo^xi) = fo^ix?) = 0. We 
accept without rigorous proof that the series converge 
and the a-dependence disappears when we approximate 
Po^i + a ) ~ ad x ip 0fl ,{xi). 

Thus, omitting the time-independent factors, we can 
express P(x,t) and F(t) in terms of eigenfunction ex- 
pansions: 

P(x 2 ,t)(xJ2^ox(x 1 )ip x(x 2 )e xt , (A2) 

A 

F(t) cx ^d^ 0AI (xi)c^ ^2)e Mt . (A3) 

A* 

Now the problem is reduced to finding a relation be- 
tween the eigenfunctions with and without absorbing 
boundaries. First, we recall that the eigenfunctions ipo\ 
without absorbing boundaries are localized in the har- 
monic region near the saddle point and contain only the 
decaying WKB "waves" , 



V>0A„ = AnC-tW/y/rfx) 

r 

4>{x) — I dzp(z) , 



p(x) = n-W(dxVo) 2 - nd^V + 2 K \. 



(A4) 
(A5) 
(A6) 



The functions are a l so localized near the sad- 
dle point. They do need to contain the growing WKB 
"waves" that cancel out the decaying "waves" at the 
boundaries x\ and x 2 . However, these additional "waves" 
decay rapidly inside the interval (x i,x 2 ) and virtually 
vanish in the saddle-point harmonic region. Therefore, 
the eigenvalues of (po^ are determined in the harmonic 
region and are the same as of the eigenfunctions ipo\ n : 
/.t = — nk. Moreover, since normalization in the har- 
monic region determines the amplitudes of the decaying 
"waves" in the WKB domain, they are equal in the cor- 
responding (fio\ n and ?Aoa„- Thus, in the WKB domain 
we obtain the following eigenfunctions with absorbing 
boundaries: 



^ 0Xn (x)=A n e-^/y^pJx) + B n e^/y^pJx), (A7) 

where A n and p{x) are the same coefficient as in the 
corresponding V'oa,, in Eq. (|A4jl . The coefficients B n 
are generally different in the two WKB sub-domains. At 
x > the amplitude B n of the growing "wave" is found 
from the boundary condition f>^\ ri (x2) = 0: 

B n e* { **> = ~A n e" t ' (x2) . (A8) 
From Eq. (|A7|) in the WKB domain we obtain 



d Mn (x) = -^pJx)(A n e-^ - B n e^ x) ) . (A9) 



where the functions tpo^ determine the eigenstates of 
the FP operator (Ce^ Vo( - x ^ K (p 0l , = lie- v °WI K ip Qli ) with 



Finally, combining Eqs. (|A4|) . (fA8|) and (|A9]) . we find 
that the WKB eigenfunctions with and without absorb- 
ing boundaries are related by 

dx<Po\ n (x 2 ) = -2p(x 2 )ip x n [x 2 ) . (A10) 
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A similar relation can be found at x±. The dependence weak-noise relation 
of p(x) on A can be neglected in the weak-noise limit, 
since the series over the eigenvalues rapidly converges 

for relevant transition times. The easiest way to verify F(t) oc P(x 2 ,t) . (All) 

that is to analyze the final expression (|30|) . Therefore, 
comparing Eqs. (|A2|1 and (|A3|) . we obtain the following 
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